"""
DaoMath Particle Physics (道数粒子物理)
从 ±0 先天均衡构造标准模型

粒子物理的 ±0 诠释:
    粒子 ↔ 反粒子 (对偶)
    费米子 ↔ 玻色子 (对偶)
    左手 ↔ 右手 (手征对偶)
    夸克 ↔ 轻子 (对偶)
    
理论层级:
    ±0 (先天均衡)
     ↓
    时空对称性 (Lorentz, Poincaré)
     ↓
    规范对称性 (U(1), SU(2), SU(3))
     ↓
    基本粒子 (夸克, 轻子, 规范玻色子)
     ↓
    相互作用 (强, 弱, 电磁, 引力)
     ↓
    Higgs 机制 (对称破缺)
     ↓
    标准模型
     ↓
    超出标准模型
"""

from __future__ import annotations
from typing import List, Tuple, Callable, Optional, Dict, Any, Union
from dataclasses import dataclass, field
from abc import ABC, abstractmethod
from enum import Enum
import numpy as np
import scipy.linalg as la
from scipy.integrate import odeint, quad
from scipy.special import gamma as gamma_function
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from core.axioms import InnateBalance, DualPair


# ============================================================================
# 物理常数
# ============================================================================

class ParticleConstants:
    """粒子物理常数"""
    C = 299792458.0                 # 光速 (m/s)
    HBAR = 1.054571817e-34          # ℏ (J·s)
    E_CHARGE = 1.602176634e-19      # 基本电荷 (C)
    M_ELECTRON = 0.51099895000      # 电子质量 (MeV/c²)
    M_PROTON = 938.27208816         # 质子质量 (MeV/c²)
    M_NEUTRON = 939.56542052        # 中子质量 (MeV/c²)
    M_W = 80379.0                   # W 玻色子质量 (MeV/c²)
    M_Z = 91187.6                   # Z 玻色子质量 (MeV/c²)
    M_HIGGS = 125100.0              # Higgs 玻色子质量 (MeV/c²)
    ALPHA_EM = 1/137.035999084      # 精细结构常数
    ALPHA_S = 0.1181                # 强耦合常数 (M_Z)
    G_F = 1.1663787e-5              # Fermi 常数 (GeV^-2)
    SIN2_THETA_W = 0.23122          # Weinberg 角
    
    # 自然单位制 (ℏ = c = 1)
    HBAR_C = HBAR * C * 1e9 / E_CHARGE  # ℏc (eV·nm)


# ============================================================================
# Level 1: 粒子分类
# ============================================================================

class Spin(Enum):
    """自旋"""
    SPIN_0 = 0
    SPIN_HALF = 0.5
    SPIN_1 = 1
    SPIN_3_HALF = 1.5
    SPIN_2 = 2


class Charge:
    """电荷"""
    def __init__(self, value: float):
        """
        Args:
            value: 电荷（以基本电荷 e 为单位）
        """
        self.value = value
    
    def __repr__(self) -> str:
        if self.value == 0:
            return "0"
        elif self.value == 1:
            return "+e"
        elif self.value == -1:
            return "-e"
        else:
            return f"{self.value:+.2f}e"
    
    def __add__(self, other: Charge) -> Charge:
        return Charge(self.value + other.value)
    
    def __neg__(self) -> Charge:
        return Charge(-self.value)


class Color(Enum):
    """色荷（强相互作用）"""
    RED = "r"
    GREEN = "g"
    BLUE = "b"
    ANTI_RED = "r̄"
    ANTI_GREEN = "ḡ"
    ANTI_BLUE = "b̄"
    COLORLESS = "colorless"  # 白色（单态）


class Generation(Enum):
    """代（世代）"""
    FIRST = 1
    SECOND = 2
    THIRD = 3


class Particle(ABC):
    """
    基本粒子抽象基类
    
    基于 ±0 的理解:
    - 每个粒子都有反粒子（正负对偶）
    - 费米子（物质）↔ 玻色子（力）
    """
    
    def __init__(self,
                 name: str,
                 mass: float,
                 charge: Charge,
                 spin: Spin):
        """
        Args:
            name: 粒子名称
            mass: 质量 (MeV/c²)
            charge: 电荷
            spin: 自旋
        """
        self.name = name
        self.mass = mass
        self.charge = charge
        self.spin = spin
        self.is_antiparticle = False
    
    @abstractmethod
    def antiparticle(self) -> Particle:
        """返回反粒子"""
        pass
    
    def is_fermion(self) -> bool:
        """判断是否为费米子（半整数自旋）"""
        return self.spin.value % 1 == 0.5
    
    def is_boson(self) -> bool:
        """判断是否为玻色子（整数自旋）"""
        return self.spin.value % 1 == 0
    
    def compton_wavelength(self) -> float:
        """
        Compton 波长
        
        λ_C = ℏ / (m c)
        
        粒子的量子尺度
        """
        if self.mass == 0:
            return float('inf')
        
        hbar_c = ParticleConstants.HBAR_C
        return hbar_c / self.mass  # nm
    
    def __repr__(self) -> str:
        anti = "̄" if self.is_antiparticle else ""
        return f"{self.name}{anti}(m={self.mass:.2f} MeV, Q={self.charge}, s={self.spin.value})"


# ============================================================================
# Level 2: 费米子（物质粒子）
# ============================================================================

class Fermion(Particle):
    """
    费米子（自旋 1/2）
    
    遵守 Pauli 不相容原理
    """
    
    def __init__(self, name: str, mass: float, charge: Charge):
        super().__init__(name, mass, charge, Spin.SPIN_HALF)


class Lepton(Fermion):
    """
    轻子
    
    不参与强相互作用
    
    体现 ±0: 带电轻子 ↔ 中微子
    """
    
    def __init__(self,
                 name: str,
                 mass: float,
                 charge: Charge,
                 generation: Generation,
                 flavor: str):
        """
        Args:
            flavor: 味（e, μ, τ, νe, νμ, ντ）
        """
        super().__init__(name, mass, charge)
        self.generation = generation
        self.flavor = flavor
    
    def antiparticle(self) -> Lepton:
        """反轻子"""
        anti = Lepton(
            name=f"anti-{self.name}",
            mass=self.mass,
            charge=-self.charge,
            generation=self.generation,
            flavor=self.flavor
        )
        anti.is_antiparticle = True
        return anti
    
    @staticmethod
    def electron() -> Lepton:
        """电子 e⁻"""
        return Lepton("electron", 0.511, Charge(-1), Generation.FIRST, "e")
    
    @staticmethod
    def muon() -> Lepton:
        """缪子 μ⁻"""
        return Lepton("muon", 105.66, Charge(-1), Generation.SECOND, "μ")
    
    @staticmethod
    def tau() -> Lepton:
        """陶子 τ⁻"""
        return Lepton("tau", 1776.86, Charge(-1), Generation.THIRD, "τ")
    
    @staticmethod
    def electron_neutrino() -> Lepton:
        """电中微子 νₑ"""
        return Lepton("electron-neutrino", 0.0, Charge(0), Generation.FIRST, "νe")
    
    @staticmethod
    def muon_neutrino() -> Lepton:
        """缪中微子 νμ"""
        return Lepton("muon-neutrino", 0.0, Charge(0), Generation.SECOND, "νμ")
    
    @staticmethod
    def tau_neutrino() -> Lepton:
        """陶中微子 ντ"""
        return Lepton("tau-neutrino", 0.0, Charge(0), Generation.THIRD, "ντ")


class Quark(Fermion):
    """
    夸克
    
    参与强相互作用，携带色荷
    
    永远禁闭在强子内（色禁闭）
    """
    
    def __init__(self,
                 name: str,
                 mass: float,
                 charge: Charge,
                 generation: Generation,
                 flavor: str,
                 color: Color = Color.RED):
        """
        Args:
            flavor: 味（u, d, c, s, t, b）
            color: 色荷
        """
        super().__init__(name, mass, charge)
        self.generation = generation
        self.flavor = flavor
        self.color = color
    
    def antiparticle(self) -> Quark:
        """反夸克"""
        # 反色荷
        color_map = {
            Color.RED: Color.ANTI_RED,
            Color.GREEN: Color.ANTI_GREEN,
            Color.BLUE: Color.ANTI_BLUE,
            Color.ANTI_RED: Color.RED,
            Color.ANTI_GREEN: Color.GREEN,
            Color.ANTI_BLUE: Color.BLUE
        }
        
        anti = Quark(
            name=f"anti-{self.name}",
            mass=self.mass,
            charge=-self.charge,
            generation=self.generation,
            flavor=self.flavor,
            color=color_map[self.color]
        )
        anti.is_antiparticle = True
        return anti
    
    @staticmethod
    def up(color: Color = Color.RED) -> Quark:
        """上夸克 u"""
        return Quark("up", 2.2, Charge(2/3), Generation.FIRST, "u", color)
    
    @staticmethod
    def down(color: Color = Color.RED) -> Quark:
        """下夸克 d"""
        return Quark("down", 4.7, Charge(-1/3), Generation.FIRST, "d", color)
    
    @staticmethod
    def charm(color: Color = Color.RED) -> Quark:
        """粲夸克 c"""
        return Quark("charm", 1270, Charge(2/3), Generation.SECOND, "c", color)
    
    @staticmethod
    def strange(color: Color = Color.RED) -> Quark:
        """奇夸克 s"""
        return Quark("strange", 95, Charge(-1/3), Generation.SECOND, "s", color)
    
    @staticmethod
    def top(color: Color = Color.RED) -> Quark:
        """顶夸克 t"""
        return Quark("top", 173000, Charge(2/3), Generation.THIRD, "t", color)
    
    @staticmethod
    def bottom(color: Color = Color.RED) -> Quark:
        """底夸克 b"""
        return Quark("bottom", 4180, Charge(-1/3), Generation.THIRD, "b", color)


# ============================================================================
# Level 3: 玻色子（力传递粒子）
# ============================================================================

class Boson(Particle):
    """
    玻色子
    
    传递相互作用
    """
    
    def __init__(self, name: str, mass: float, charge: Charge, spin: Spin):
        super().__init__(name, mass, charge, spin)


class GaugeBoson(Boson):
    """
    规范玻色子（自旋 1）
    
    传递基本相互作用
    
    基于 ±0 的理解:
    - 规范对称性 → 守恒定律
    - 局域对称性 → 相互作用
    """
    
    def __init__(self, 
                 name: str,
                 mass: float,
                 charge: Charge,
                 interaction: str):
        """
        Args:
            interaction: 相互作用类型（electromagnetic, weak, strong）
        """
        super().__init__(name, mass, charge, Spin.SPIN_1)
        self.interaction = interaction
    
    def antiparticle(self) -> GaugeBoson:
        """规范玻色子的反粒子"""
        anti = GaugeBoson(
            name=f"anti-{self.name}",
            mass=self.mass,
            charge=-self.charge,
            interaction=self.interaction
        )
        anti.is_antiparticle = True
        return anti
    
    @staticmethod
    def photon() -> GaugeBoson:
        """
        光子 γ
        
        电磁相互作用
        无质量，电中性
        """
        return GaugeBoson("photon", 0.0, Charge(0), "electromagnetic")
    
    @staticmethod
    def gluon(color_charge: str = "rḡ") -> GaugeBoson:
        """
        胶子 g
        
        强相互作用
        携带色荷（色-反色组合）
        8 种独立态
        """
        g = GaugeBoson("gluon", 0.0, Charge(0), "strong")
        g.color_charge = color_charge
        return g
    
    @staticmethod
    def w_plus() -> GaugeBoson:
        """W⁺ 玻色子（弱相互作用）"""
        return GaugeBoson("W+", ParticleConstants.M_W, Charge(1), "weak")
    
    @staticmethod
    def w_minus() -> GaugeBoson:
        """W⁻ 玻色子（弱相互作用）"""
        return GaugeBoson("W-", ParticleConstants.M_W, Charge(-1), "weak")
    
    @staticmethod
    def z_boson() -> GaugeBoson:
        """Z⁰ 玻色子（弱相互作用）"""
        return GaugeBoson("Z0", ParticleConstants.M_Z, Charge(0), "weak")


class HiggsBoson(Boson):
    """
    Higgs 玻色子（自旋 0）
    
    赋予粒子质量
    
    基于 ±0 的理解:
    - Higgs 场真空期望值 ≠ 0
    - 自发对称破缺
    - 规范对称 ↔ 质量生成
    """
    
    def __init__(self):
        super().__init__(
            "Higgs",
            ParticleConstants.M_HIGGS,
            Charge(0),
            Spin.SPIN_0
        )
        self.vev = 246000.0  # 真空期望值 (MeV)
    
    def antiparticle(self) -> HiggsBoson:
        """Higgs 玻色子是自己的反粒子"""
        return self
    
    def yukawa_coupling(self, fermion: Fermion) -> float:
        """
        Yukawa 耦合常数
        
        y_f = √2 m_f / v
        
        决定费米子质量
        """
        return np.sqrt(2) * fermion.mass / self.vev


# ============================================================================
# Level 4: 复合粒子（强子）
# ============================================================================

class Hadron(Particle):
    """
    强子
    
    由夸克组成的复合粒子
    
    - 介子：夸克-反夸克对（q q̄）
    - 重子：三夸克（qqq）
    """
    
    def __init__(self,
                 name: str,
                 mass: float,
                 charge: Charge,
                 spin: Spin,
                 quark_content: List[Quark]):
        super().__init__(name, mass, charge, spin)
        self.quarks = quark_content
        self.baryon_number = self._calculate_baryon_number()
    
    def _calculate_baryon_number(self) -> int:
        """
        重子数
        
        夸克: +1/3
        反夸克: -1/3
        """
        B = 0
        for q in self.quarks:
            B += 1 if not q.is_antiparticle else -1
        return B // 3
    
    def is_meson(self) -> bool:
        """判断是否为介子"""
        return self.baryon_number == 0
    
    def is_baryon(self) -> bool:
        """判断是否为重子"""
        return abs(self.baryon_number) == 1
    
    def antiparticle(self) -> Hadron:
        """反强子"""
        anti_quarks = [q.antiparticle() for q in self.quarks]
        
        anti = Hadron(
            name=f"anti-{self.name}",
            mass=self.mass,
            charge=-self.charge,
            spin=self.spin,
            quark_content=anti_quarks
        )
        anti.is_antiparticle = True
        return anti


class Meson(Hadron):
    """
    介子（q q̄）
    
    整数自旋（玻色子）
    """
    
    def __init__(self,
                 name: str,
                 mass: float,
                 charge: Charge,
                 spin: Spin,
                 quark: Quark,
                 antiquark: Quark):
        super().__init__(name, mass, charge, spin, [quark, antiquark])
    
    @staticmethod
    def pion_plus() -> Meson:
        """
        π⁺ 介子
        
        u d̄
        """
        u = Quark.up()
        d_bar = Quark.down().antiparticle()
        
        return Meson("π+", 139.57, Charge(1), Spin.SPIN_0, u, d_bar)
    
    @staticmethod
    def pion_minus() -> Meson:
        """
        π⁻ 介子
        
        d ū
        """
        d = Quark.down()
        u_bar = Quark.up().antiparticle()
        
        return Meson("π-", 139.57, Charge(-1), Spin.SPIN_0, d, u_bar)
    
    @staticmethod
    def pion_zero() -> Meson:
        """
        π⁰ 介子
        
        (u ū - d d̄) / √2
        """
        u = Quark.up()
        u_bar = u.antiparticle()
        
        return Meson("π0", 134.98, Charge(0), Spin.SPIN_0, u, u_bar)
    
    @staticmethod
    def kaon_plus() -> Meson:
        """
        K⁺ 介子
        
        u s̄
        """
        u = Quark.up()
        s_bar = Quark.strange().antiparticle()
        
        return Meson("K+", 493.68, Charge(1), Spin.SPIN_0, u, s_bar)


class Baryon(Hadron):
    """
    重子（qqq）
    
    半整数自旋（费米子）
    """
    
    def __init__(self,
                 name: str,
                 mass: float,
                 charge: Charge,
                 spin: Spin,
                 quarks: List[Quark]):
        if len(quarks) != 3:
            raise ValueError("Baryon must have exactly 3 quarks")
        
        super().__init__(name, mass, charge, spin, quarks)
    
    @staticmethod
    def proton() -> Baryon:
        """
        质子 p
        
        uud
        """
        u1 = Quark.up(Color.RED)
        u2 = Quark.up(Color.GREEN)
        d = Quark.down(Color.BLUE)
        
        return Baryon("proton", 938.27, Charge(1), Spin.SPIN_HALF, [u1, u2, d])
    
    @staticmethod
    def neutron() -> Baryon:
        """
        中子 n
        
        udd
        """
        u = Quark.up(Color.RED)
        d1 = Quark.down(Color.GREEN)
        d2 = Quark.down(Color.BLUE)
        
        return Baryon("neutron", 939.57, Charge(0), Spin.SPIN_HALF, [u, d1, d2])
    
    @staticmethod
    def lambda_baryon() -> Baryon:
        """
        Λ 重子
        
        uds
        """
        u = Quark.up(Color.RED)
        d = Quark.down(Color.GREEN)
        s = Quark.strange(Color.BLUE)
        
        return Baryon("Λ", 1115.68, Charge(0), Spin.SPIN_HALF, [u, d, s])


# ============================================================================
# Level 5: 相互作用
# ============================================================================

class Interaction(ABC):
    """
    基本相互作用抽象基类
    
    基于 ±0 的理解:
    - 吸引 ↔ 排斥
    - 虚粒子交换：正反对偶
    """
    
    def __init__(self, 
                 name: str,
                 mediator: GaugeBoson,
                 coupling_constant: float,
                 range: float):
        """
        Args:
            name: 相互作用名称
            mediator: 传递粒子
            coupling_constant: 耦合常数
            range: 作用范围 (m)
        """
        self.name = name
        self.mediator = mediator
        self.alpha = coupling_constant
        self.range = range
    
    @abstractmethod
    def potential(self, r: float) -> float:
        """
        势能 V(r)
        
        Args:
            r: 距离 (m)
        """
        pass
    
    @abstractmethod
    def cross_section(self, energy: float) -> float:
        """
        散射截面 σ(E)
        
        Args:
            energy: 能量 (MeV)
        """
        pass


class ElectromagneticInteraction(Interaction):
    """
    电磁相互作用
    
    - 传递粒子：光子 γ
    - 规范对称性：U(1)
    - 长程力
    
    基于 ±0 的理解:
    - 同号电荷排斥（正正/负负）
    - 异号电荷吸引（正负）
    """
    
    def __init__(self):
        super().__init__(
            name="electromagnetic",
            mediator=GaugeBoson.photon(),
            coupling_constant=ParticleConstants.ALPHA_EM,
            range=float('inf')  # 长程
        )
    
    def potential(self, r: float) -> float:
        """
        Coulomb 势
        
        V(r) = α ℏc / r
        
        α = e²/(4πε₀ℏc) ≈ 1/137
        """
        hbar_c = ParticleConstants.HBAR_C
        
        return self.alpha * hbar_c / r  # MeV·nm / nm = MeV
    
    def cross_section(self, energy: float) -> float:
        """
        Rutherford 散射截面
        
        dσ/dΩ = (α ℏc / 2E)² / sin⁴(θ/2)
        """
        # 简化：总截面（积分）
        hbar_c = ParticleConstants.HBAR_C
        
        sigma_0 = (self.alpha * hbar_c / energy)**2
        
        return sigma_0  # 相对单位
    
    def fine_structure_splitting(self, n: int) -> float:
        """
        精细结构分裂
        
        ΔE ~ α² m_e c² / n³
        
        体现相对论修正
        """
        m_e = ParticleConstants.M_ELECTRON
        
        return self.alpha**2 * m_e / n**3


class WeakInteraction(Interaction):
    """
    弱相互作用
    
    - 传递粒子：W±, Z⁰
    - 规范对称性：SU(2)_L × U(1)_Y
    - 短程力
    - 破坏宇称（P）、电荷共轭（C）
    
    基于 ±0 的理解:
    - 左手 ↔ 右手（手征不对称）
    - β 衰变：中子 → 质子（味改变）
    """
    
    def __init__(self):
        # 使用 W 玻色子作为代表
        super().__init__(
            name="weak",
            mediator=GaugeBoson.w_plus(),
            coupling_constant=ParticleConstants.G_F,  # Fermi 常数
            range=2e-18  # ~ 0.002 fm
        )
        self.M_W = ParticleConstants.M_W
        self.M_Z = ParticleConstants.M_Z
    
    def potential(self, r: float) -> float:
        """
        Yukawa 势
        
        V(r) = -g² ℏc e^(-m_W c r / ℏ) / (4π r)
        
        有质量传递粒子 → 短程
        """
        hbar_c = ParticleConstants.HBAR_C
        
        # 转换单位
        m_W = self.M_W / hbar_c  # nm^-1
        
        return -self.alpha * hbar_c * np.exp(-m_W * r) / r
    
    def cross_section(self, energy: float) -> float:
        """
        弱相互作用散射截面
        
        σ ~ G_F² s
        
        其中 s = (E_CM)² 是质心系能量平方
        """
        # Fermi 常数（GeV^-2）
        G_F = self.alpha
        
        # 能量（MeV → GeV）
        E_GeV = energy / 1000
        
        sigma = G_F**2 * E_GeV**2
        
        return sigma  # GeV^-2 ~ 0.4 mb
    
    def beta_decay_rate(self, Q_value: float) -> float:
        """
        β 衰变率
        
        Γ ~ G_F² Q⁵
        
        Args:
            Q_value: 衰变能量 (MeV)
        """
        G_F = self.alpha
        
        return G_F**2 * Q_value**5  # s^-1
    
    def neutrino_oscillation_probability(self,
                                        L: float,
                                        E: float,
                                        delta_m_sq: float,
                                        theta: float) -> float:
        """
        中微子振荡概率
        
        P(ν_α → ν_β) = sin²(2θ) sin²(1.27 Δm² L / E)
        
        体现 ±0: 味本征态 ↔ 质量本征态的混合
        
        Args:
            L: 传播距离 (km)
            E: 能量 (GeV)
            delta_m_sq: 质量平方差 (eV²)
            theta: 混合角
        """
        # 振荡相位
        phase = 1.27 * delta_m_sq * L / E
        
        P = np.sin(2 * theta)**2 * np.sin(phase)**2
        
        return P


class StrongInteraction(Interaction):
    """
    强相互作用（量子色动力学 QCD）
    
    - 传递粒子：胶子 g
    - 规范对称性：SU(3)_color
    - 短程力（色禁闭）
    - 渐近自由
    
    基于 ±0 的理解:
    - 色荷：红-反红，绿-反绿，蓝-反蓝
    - 夸克 ↔ 反夸克
    - 胶子自相互作用
    """
    
    def __init__(self):
        super().__init__(
            name="strong",
            mediator=GaugeBoson.gluon(),
            coupling_constant=ParticleConstants.ALPHA_S,
            range=1e-15  # ~ 1 fm
        )
    
    def potential(self, r: float) -> float:
        """
        Cornell 势（唯象）
        
        V(r) = -4α_s/(3r) + kr
        
        - 短程：Coulomb 型（渐近自由）
        - 长程：线性禁闭
        
        k ~ 1 GeV/fm 是弦张力
        """
        hbar_c = ParticleConstants.HBAR_C
        k = 1000.0  # MeV/fm
        
        r_fm = r * 1e6  # nm → fm
        
        V_coulomb = -4 * self.alpha * hbar_c / (3 * r_fm)
        V_linear = k * r_fm
        
        return V_coulomb + V_linear
    
    def cross_section(self, energy: float) -> float:
        """
        强相互作用截面
        
        σ ~ (ℏc / E)² × α_s²
        """
        hbar_c = ParticleConstants.HBAR_C
        
        sigma = (hbar_c / energy)**2 * self.alpha**2
        
        return sigma
    
    def running_coupling(self, Q: float) -> float:
        """
        流动耦合常数（渐近自由）
        
        α_s(Q²) = α_s(M_Z²) / [1 + (33-2n_f)α_s(M_Z²)/(12π) ln(Q²/M_Z²)]
        
        Q → ∞: α_s → 0 (渐近自由)
        Q → Λ_QCD: α_s → ∞ (色禁闭)
        
        Args:
            Q: 能量标度 (MeV)
        """
        M_Z = ParticleConstants.M_Z
        alpha_s_MZ = self.alpha
        
        n_f = 5  # 有效夸克数（忽略顶夸克）
        
        beta0 = (33 - 2 * n_f) / (12 * np.pi)
        
        log_term = np.log(Q**2 / M_Z**2)
        
        alpha_s_Q = alpha_s_MZ / (1 + beta0 * alpha_s_MZ * log_term)
        
        return alpha_s_Q
    
    def confinement_scale(self) -> float:
        """
        色禁闭标度 Λ_QCD
        
        ~ 200 MeV
        """
        return 200.0  # MeV


# ============================================================================
# Level 6: 标准模型
# ============================================================================

class StandardModel:
    """
    标准模型
    
    统一描述基本粒子和相互作用（除引力外）
    
    规范群：SU(3)_C × SU(2)_L × U(1)_Y
    
    基于 ±0 的理解:
    - 粒子-反粒子对称（C）
    - 左-右手对称性破缺（弱相互作用）
    - Higgs 机制：对称破缺 → 质量生成
    """
    
    def __init__(self):
        # 费米子（3 代）
        self.leptons = {
            1: [Lepton.electron(), Lepton.electron_neutrino()],
            2: [Lepton.muon(), Lepton.muon_neutrino()],
            3: [Lepton.tau(), Lepton.tau_neutrino()]
        }
        
        self.quarks = {
            1: [Quark.up(), Quark.down()],
            2: [Quark.charm(), Quark.strange()],
            3: [Quark.top(), Quark.bottom()]
        }
        
        # 规范玻色子
        self.gauge_bosons = {
            'photon': GaugeBoson.photon(),
            'W+': GaugeBoson.w_plus(),
            'W-': GaugeBoson.w_minus(),
            'Z0': GaugeBoson.z_boson(),
            'gluon': [GaugeBoson.gluon(f"g{i}") for i in range(8)]
        }
        
        # Higgs 玻色子
        self.higgs = HiggsBoson()
        
        # 相互作用
        self.interactions = {
            'electromagnetic': ElectromagneticInteraction(),
            'weak': WeakInteraction(),
            'strong': StrongInteraction()
        }
    
    def particle_count(self) -> Dict[str, int]:
        """
        标准模型粒子计数
        
        Returns:
            粒子数统计
        """
        count = {
            'quarks': 6,           # u, d, c, s, t, b
            'leptons': 6,          # e, μ, τ, νe, νμ, ντ
            'gauge_bosons': 12,    # γ, W±, Z, 8g
            'higgs': 1,
            'total_fermions': 12,
            'total_bosons': 13,
            'total': 25
        }
        
        return count
    
    def symmetry_groups(self) -> Dict[str, str]:
        """
        规范对称性
        
        SU(3)_C: 色 SU(3)（强相互作用）
        SU(2)_L: 弱同位旋 SU(2)（左手费米子）
        U(1)_Y: 超荷 U(1)（弱超荷）
        """
        return {
            'strong': 'SU(3)_C',
            'electroweak': 'SU(2)_L × U(1)_Y',
            'full': 'SU(3)_C × SU(2)_L × U(1)_Y'
        }
    
    def electroweak_unification(self) -> Dict[str, Any]:
        """
        电弱统一
        
        Glashow-Weinberg-Salam 理论
        
        电磁 + 弱 → 电弱
        
        对称破缺后：
        SU(2)_L × U(1)_Y → U(1)_EM
        """
        sin2_theta_W = ParticleConstants.SIN2_THETA_W
        cos2_theta_W = 1 - sin2_theta_W
        
        return {
            'unified_symmetry': 'SU(2)_L × U(1)_Y',
            'broken_symmetry': 'U(1)_EM',
            'Weinberg_angle': np.arcsin(np.sqrt(sin2_theta_W)),
            'W_mass': ParticleConstants.M_W,
            'Z_mass': ParticleConstants.M_Z,
            'mass_ratio': ParticleConstants.M_W / ParticleConstants.M_Z,
            'prediction': np.sqrt(cos2_theta_W)
        }
    
    def higgs_mechanism(self) -> Dict[str, Any]:
        """
        Higgs 机制
        
        自发对称破缺 → 质量生成
        
        体现 ±0: 真空期望值 ⟨φ⟩ ≠ 0 打破对称
        """
        v = self.higgs.vev
        
        # W, Z 质量关系
        g2 = 0.653  # SU(2) 耦合
        g1 = 0.357  # U(1) 耦合
        
        M_W_theory = g2 * v / 2
        M_Z_theory = np.sqrt(g2**2 + g1**2) * v / 2
        
        return {
            'vacuum_expectation': v,
            'W_mass_theory': M_W_theory,
            'W_mass_exp': ParticleConstants.M_W,
            'Z_mass_theory': M_Z_theory,
            'Z_mass_exp': ParticleConstants.M_Z,
            'goldstone_bosons': 3,  # 被 W±, Z 吃掉
            'higgs_remaining': 1
        }
    
    def fermion_masses(self) -> Dict[str, Dict[str, float]]:
        """
        费米子质量（通过 Yukawa 耦合）
        
        m_f = y_f v / √2
        """
        masses = {
            'quarks': {},
            'leptons': {}
        }
        
        for gen in [1, 2, 3]:
            for q in self.quarks[gen]:
                masses['quarks'][q.name] = q.mass
            
            for l in self.leptons[gen]:
                masses['leptons'][l.name] = l.mass
        
        return masses
    
    def ckm_matrix(self) -> np.ndarray:
        """
        CKM 矩阵（夸克混合）
        
        描述不同代夸克间的弱相互作用
        
        体现 ±0: 味本征态 ↔ 质量本征态
        """
        # 简化的 CKM 矩阵（实验值）
        V_CKM = np.array([
            [0.97435, 0.22500, 0.00369],
            [0.22486, 0.97349, 0.04182],
            [0.00857, 0.04110, 0.99915]
        ])
        
        # 验证幺正性
        unitarity_check = V_CKM @ V_CKM.T.conj()
        
        return V_CKM
    
    def pmns_matrix(self) -> np.ndarray:
        """
        PMNS 矩阵（轻子混合）
        
        描述中微子振荡
        """
        # 简化的 PMNS 矩阵
        theta12 = 33.82 * np.pi / 180
        theta23 = 48.3 * np.pi / 180
        theta13 = 8.61 * np.pi / 180
        
        s12, c12 = np.sin(theta12), np.cos(theta12)
        s23, c23 = np.sin(theta23), np.cos(theta23)
        s13, c13 = np.sin(theta13), np.cos(theta13)
        
        U_PMNS = np.array([
            [c12*c13, s12*c13, s13],
            [-s12*c23 - c12*s23*s13, c12*c23 - s12*s23*s13, s23*c13],
            [s12*s23 - c12*c23*s13, -c12*s23 - s12*c23*s13, c23*c13]
        ])
        
        return U_PMNS


class FeynmanDiagram:
    """
    Feynman 图
    
    描述粒子相互作用的图形表示
    
    基于 ±0 的理解:
    - 顶点：相互作用
    - 传播子：虚粒子（正反对）
    - 外线：实粒子
    """
    
    def __init__(self, 
                 process: str,
                 initial_particles: List[Particle],
                 final_particles: List[Particle],
                 mediator: Optional[GaugeBoson] = None):
        """
        Args:
            process: 过程描述
            initial_particles: 初态粒子
            final_particles: 末态粒子
            mediator: 传递粒子
        """
        self.process = process
        self.initial = initial_particles
        self.final = final_particles
        self.mediator = mediator
    
    def conservation_check(self) -> Dict[str, bool]:
        """
        守恒律检验
        
        - 能量-动量守恒
        - 电荷守恒
        - 重子数守恒
        - 轻子数守恒
        """
        checks = {}
        
        # 电荷守恒
        Q_initial = sum(p.charge.value for p in self.initial)
        Q_final = sum(p.charge.value for p in self.final)
        checks['charge'] = abs(Q_initial - Q_final) < 1e-10
        
        # 重子数守恒（如果适用）
        B_initial = sum(getattr(p, 'baryon_number', 0) for p in self.initial)
        B_final = sum(getattr(p, 'baryon_number', 0) for p in self.final)
        checks['baryon_number'] = B_initial == B_final
        
        return checks
    
    @staticmethod
    def compton_scattering() -> FeynmanDiagram:
        """
        Compton 散射
        
        γ + e⁻ → γ + e⁻
        """
        photon = GaugeBoson.photon()
        electron = Lepton.electron()
        
        return FeynmanDiagram(
            process="Compton scattering",
            initial_particles=[photon, electron],
            final_particles=[photon, electron],
            mediator=electron
        )
    
    @staticmethod
    def beta_decay() -> FeynmanDiagram:
        """
        β⁻ 衰变
        
        n → p + e⁻ + ν̄ₑ
        
        夸克层面: d → u + W⁻ → u + e⁻ + ν̄ₑ
        """
        neutron = Baryon.neutron()
        proton = Baryon.proton()
        electron = Lepton.electron()
        neutrino = Lepton.electron_neutrino().antiparticle()
        
        return FeynmanDiagram(
            process="β⁻ decay",
            initial_particles=[neutron],
            final_particles=[proton, electron, neutrino],
            mediator=GaugeBoson.w_minus()
        )
    
    @staticmethod
    def pair_production() -> FeynmanDiagram:
        """
        正负电子对产生
        
        γ → e⁺ + e⁻
        """
        photon = GaugeBoson.photon()
        electron = Lepton.electron()
        positron = electron.antiparticle()
        
        return FeynmanDiagram(
            process="pair production",
            initial_particles=[photon],
            final_particles=[electron, positron]
        )


# ============================================================================
# Level 7: 超出标准模型
# ============================================================================

class BeyondStandardModel:
    """
    超出标准模型（BSM）
    
    标准模型的未解问题：
    - 暗物质
    - 中微子质量
    - CP 破坏
    - 层次问题
    - 引力量子化
    """
    
    @staticmethod
    def dark_matter_candidates() -> List[str]:
        """
        暗物质候选
        
        体现 ±0: 可见物质 ↔ 暗物质
        """
        return [
            "WIMPs (Weakly Interacting Massive Particles)",
            "Axions",
            "Sterile neutrinos",
            "Primordial black holes"
        ]
    
    @staticmethod
    def neutrino_mass_mechanisms() -> List[str]:
        """
        中微子质量机制
        
        标准模型预言中微子无质量，但实验观测到振荡
        """
        return [
            "Seesaw mechanism",
            "Radiative mass generation",
            "Dirac mass (right-handed neutrinos)"
        ]
    
    @staticmethod
    def grand_unification() -> Dict[str, Any]:
        """
        大统一理论（GUT）
        
        SU(3) × SU(2) × U(1) → SU(5) / SO(10) / E₆
        
        统一所有相互作用（除引力）
        """
        return {
            'unification_scale': 1e16,  # GeV
            'proton_decay': 'p → e⁺ + π⁰',
            'lifetime': '>1e34 years',
            'gauge_groups': ['SU(5)', 'SO(10)', 'E₆']
        }
    
    @staticmethod
    def supersymmetry() -> Dict[str, str]:
        """
        超对称（SUSY）
        
        费米子 ↔ 玻色子
        
        体现 ±0: 每个粒子都有超对称伙伴
        """
        partners = {
            'electron': 'selectron (ẽ)',
            'quark': 'squark (q̃)',
            'photon': 'photino (γ̃)',
            'gluon': 'gluino (g̃)',
            'Higgs': 'Higgsino (H̃)',
            'W/Z': 'Wino/Zino (W̃/Z̃)'
        }
        
        return partners
    
    @staticmethod
    def extra_dimensions() -> Dict[str, Any]:
        """
        额外维度
        
        解决层次问题
        
        Kaluza-Klein 理论
        """
        return {
            'compactification_scale': 'TeV - Planck',
            'KK_modes': 'tower of massive particles',
            'gravity_weakness': 'diluted in extra dimensions'
        }


# ============================================================================
# 综合演示系统
# ============================================================================

class ParticlePhysicsDemo:
    """粒子物理综合演示"""
    
    @staticmethod
    def demonstrate_particles():
        """演示基本粒子"""
        print("\n" + "=" * 70)
        print("标准模型基本粒子")
        print("=" * 70)
        
        print("\n【轻子】")
        print("第一代:")
        e = Lepton.electron()
        nu_e = Lepton.electron_neutrino()
        print(f"  {e}")
        print(f"  {nu_e}")
        
        print("\n第二代:")
        mu = Lepton.muon()
        nu_mu = Lepton.muon_neutrino()
        print(f"  {mu}")
        print(f"  {nu_mu}")
        
        print("\n第三代:")
        tau = Lepton.tau()
        nu_tau = Lepton.tau_neutrino()
        print(f"  {tau}")
        print(f"  {nu_tau}")
        
        print("\n【夸克】")
        print("第一代:")
        u = Quark.up()
        d = Quark.down()
        print(f"  {u}")
        print(f"  {d}")
        
        print("\n第二代:")
        c = Quark.charm()
        s = Quark.strange()
        print(f"  {c}")
        print(f"  {s}")
        
        print("\n第三代:")
        t = Quark.top()
        b = Quark.bottom()
        print(f"  {t}")
        print(f"  {b}")
        
        print("\n【规范玻色子】")
        photon = GaugeBoson.photon()
        W_plus = GaugeBoson.w_plus()
        Z = GaugeBoson.z_boson()
        gluon = GaugeBoson.gluon()
        
        print(f"  {photon} (电磁)")
        print(f"  {W_plus} (弱)")
        print(f"  {Z} (弱)")
        print(f"  {gluon} ×8 (强)")
        
        print("\n【Higgs 玻色子】")
        higgs = HiggsBoson()
        print(f"  {higgs}")
        print(f"  真空期望值 v = {higgs.vev:.0f} MeV")
    
    @staticmethod
    def demonstrate_hadrons():
        """演示强子"""
        print("\n" + "=" * 70)
        print("强子（夸克复合粒子）")
        print("=" * 70)
        
        print("\n【介子 (q q̄)】")
        pi_plus = Meson.pion_plus()
        pi_minus = Meson.pion_minus()
        pi_zero = Meson.pion_zero()
        K_plus = Meson.kaon_plus()
        
        print(f"  {pi_plus}")
        print(f"  {pi_minus}")
        print(f"  {pi_zero}")
        print(f"  {K_plus}")
        
        print("\n【重子 (qqq)】")
        proton = Baryon.proton()
        neutron = Baryon.neutron()
        lambda_b = Baryon.lambda_baryon()
        
        print(f"  {proton}")
        print(f"    夸克组成: {[q.name for q in proton.quarks]}")
        print(f"    重子数: {proton.baryon_number}")
        
        print(f"  {neutron}")
        print(f"    夸克组成: {[q.name for q in neutron.quarks]}")
        print(f"    重子数: {neutron.baryon_number}")
        
        print(f"  {lambda_b}")
        print(f"    夸克组成: {[q.name for q in lambda_b.quarks]}")
    
    @staticmethod
    def demonstrate_interactions():
        """演示相互作用"""
        print("\n" + "=" * 70)
        print("基本相互作用")
        print("=" * 70)
        
        # 电磁相互作用
        em = ElectromagneticInteraction()
        print(f"\n【电磁相互作用】")
        print(f"  传递粒子: {em.mediator.name}")
        print(f"  耦合常数: α = {em.alpha:.6f}")
        print(f"  作用范围: ∞")
        
        r_vals = [0.1, 1.0, 10.0]  # nm
        print(f"\n  Coulomb 势:")
        for r in r_vals:
            V = em.potential(r)
            print(f"    r = {r:.1f} nm: V = {V:.4f} MeV")
        
        # 弱相互作用
        weak = WeakInteraction()
        print(f"\n【弱相互作用】")
        print(f"  传递粒子: W±, Z⁰")
        print(f"  Fermi 常数: G_F = {weak.alpha:.4e} GeV^-2")
        print(f"  作用范围: {weak.range:.2e} m")
        
        print(f"\n  中微子振荡:")
        L = 295  # km (Daya Bay)
        E = 3.0  # MeV
        delta_m_sq = 2.5e-3  # eV²
        theta = 33.82 * np.pi / 180
        
        P = weak.neutrino_oscillation_probability(L, E/1000, delta_m_sq, theta)
        print(f"    P(νₑ → νₑ) = {P:.4f}")
        
        # 强相互作用
        strong = StrongInteraction()
        print(f"\n【强相互作用】")
        print(f"  传递粒子: 8 种胶子")
        print(f"  耦合常数: α_s(M_Z) = {strong.alpha:.4f}")
        print(f"  作用范围: {strong.range:.2e} m")
        
        print(f"\n  渐近自由:")
        Q_vals = [1000, 10000, 91187, 173000]  # MeV
        for Q in Q_vals:
            alpha_s_Q = strong.running_coupling(Q)
            print(f"    Q = {Q/1000:.1f} GeV: α_s = {alpha_s_Q:.4f}")
        
        print(f"\n  色禁闭标度:")
        Lambda_QCD = strong.confinement_scale()
        print(f"    Λ_QCD ≈ {Lambda_QCD:.0f} MeV")
    
    @staticmethod
    def demonstrate_standard_model():
        """演示标准模型"""
        print("\n" + "=" * 70)
        print("标准模型")
        print("=" * 70)
        
        sm = StandardModel()
        
        # 粒子计数
        count = sm.particle_count()
        print(f"\n粒子计数:")
        print(f"  夸克: {count['quarks']} 种 × 3 色 = 18")
        print(f"  轻子: {count['leptons']} 种")
        print(f"  规范玻色子: {count['gauge_bosons']} 种")
        print(f"  Higgs: {count['higgs']} 种")
        print(f"  总计: {count['total']} 种基本粒子")
        
        # 规范对称性
        symmetries = sm.symmetry_groups()
        print(f"\n规范对称性:")
        print(f"  强相互作用: {symmetries['strong']}")
        print(f"  电弱统一: {symmetries['electroweak']}")
        print(f"  完整群: {symmetries['full']}")
        
        # 电弱统一
        ew = sm.electroweak_unification()
        print(f"\n电弱统一:")
        print(f"  Weinberg 角: θ_W = {ew['Weinberg_angle']*180/np.pi:.2f}°")
        print(f"  sin²θ_W = {ParticleConstants.SIN2_THETA_W:.5f}")
        print(f"  M_W/M_Z = {ew['mass_ratio']:.5f}")
        print(f"  理论预言: cos θ_W = {ew['prediction']:.5f}")
        
        # Higgs 机制
        higgs_info = sm.higgs_mechanism()
        print(f"\nHiggs 机制:")
        print(f"  真空期望值: v = {higgs_info['vacuum_expectation']:.0f} MeV")
        print(f"  W 质量（理论）: {higgs_info['W_mass_theory']:.0f} MeV")
        print(f"  W 质量（实验）: {higgs_info['W_mass_exp']:.0f} MeV")
        print(f"  Z 质量（理论）: {higgs_info['Z_mass_theory']:.0f} MeV")
        print(f"  Z 质量（实验）: {higgs_info['Z_mass_exp']:.0f} MeV")
        
        # CKM 矩阵
        V_CKM = sm.ckm_matrix()
        print(f"\nCKM 矩阵（夸克混合）:")
        print(V_CKM)
        
        # PMNS 矩阵
        U_PMNS = sm.pmns_matrix()
        print(f"\nPMNS 矩阵（轻子混合）:")
        print(U_PMNS)
    
    @staticmethod
    def demonstrate_feynman_diagrams():
        """演示 Feynman 图"""
        print("\n" + "=" * 70)
        print("Feynman 图与粒子过程")
        print("=" * 70)
        
        # Compton 散射
        compton = FeynmanDiagram.compton_scattering()
        print(f"\n{compton.process}:")
        print(f"  初态: {[p.name for p in compton.initial]}")
        print(f"  末态: {[p.name for p in compton.final]}")
        print(f"  守恒律: {compton.conservation_check()}")
        
        # β 衰变
        beta = FeynmanDiagram.beta_decay()
        print(f"\n{beta.process}:")
        print(f"  初态: {[p.name for p in beta.initial]}")
        print(f"  末态: {[p.name for p in beta.final]}")
        print(f"  传递粒子: {beta.mediator.name}")
        print(f"  守恒律: {beta.conservation_check()}")
        
        # 正负电子对产生
        pair = FeynmanDiagram.pair_production()
        print(f"\n{pair.process}:")
        print(f"  初态: {[p.name for p in pair.initial]}")
        print(f"  末态: {[p.name for p in pair.final]}")
        print(f"  阈值能量: 2m_e c² = {2*ParticleConstants.M_ELECTRON:.3f} MeV")
    
    @staticmethod
    def demonstrate_bsm():
        """演示超出标准模型"""
        print("\n" + "=" * 70)
        print("超出标准模型 (BSM)")
        print("=" * 70)
        
        # 暗物质
        print("\n暗物质候选:")
        for candidate in BeyondStandardModel.dark_matter_candidates():
            print(f"  • {candidate}")
        
        # 中微子质量
        print("\n中微子质量机制:")
        for mechanism in BeyondStandardModel.neutrino_mass_mechanisms():
            print(f"  • {mechanism}")
        
        # 大统一
        gut = BeyondStandardModel.grand_unification()
        print(f"\n大统一理论 (GUT):")
        print(f"  统一标度: M_GUT ~ {gut['unification_scale']:.0e} GeV")
        print(f"  质子衰变: {gut['proton_decay']}")
        print(f"  质子寿命: τ_p > {gut['lifetime']}")
        
        # 超对称
        print(f"\n超对称伙伴:")
        susy = BeyondStandardModel.supersymmetry()
        for particle, partner in list(susy.items())[:3]:
            print(f"  {particle} ↔ {partner}")
        
        # 额外维度
        extra_dim = BeyondStandardModel.extra_dimensions()
        print(f"\n额外维度:")
        print(f"  紧致化标度: {extra_dim['compactification_scale']}")
        print(f"  KK 模式: {extra_dim['KK_modes']}")
        print(f"  引力弱化: {extra_dim['gravity_weakness']}")


def main():
    """主函数：完整演示"""
    
    print("=" * 70)
    print("DaoMath Particle Physics")
    print("从 ±0 构造标准模型")
    print("=" * 70)
    
    # 1. 基本粒子
    ParticlePhysicsDemo.demonstrate_particles()
    
    # 2. 强子
    ParticlePhysicsDemo.demonstrate_hadrons()
    
    # 3. 相互作用
    ParticlePhysicsDemo.demonstrate_interactions()
    
    # 4. 标准模型
    ParticlePhysicsDemo.demonstrate_standard_model()
    
    # 5. Feynman 图
    ParticlePhysicsDemo.demonstrate_feynman_diagrams()
    
    # 6. 超出标准模型
    ParticlePhysicsDemo.demonstrate_bsm()
    
    print("\n" + "=" * 70)
    print("粒子物理构造完成！")
    print("=" * 70)


if __name__ == "__main__":
    main()
